import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint

# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei']  # 使用 SimHei 字体
plt.rcParams['axes.unicode_minus'] = False  # 解决负号显示问题

#命名规则
#Natural growth rate of crops --> ngrc
#
# 定义微分方程组
def model(y, t, params):
    C, W, I, B, Bt = y  #设置初始种群数量
    r_C, a_CI, d_C, r_W, h, d_W, e_IC, a_IC, m_IB, m_IBt, p, e_BI, a_IB, d_B, e_BtI, d_Bt = params

    #r_C:农作物的自然增长率            #a_CI:昆虫对农作物的攻击率
    #d_C:自然死亡率                   #r_W:杂草的自然增长率
    #h:除草剂的效果系数                #d_W:自然死亡率
    #e_IC:昆虫从农作物获得能量的效率    #a_IC:昆虫捕食农作物的效率
    #m_IB:鸟类捕食昆虫的效率           #p:杀虫剂的效果系数
    #e_BI:鸟类从昆虫获得能量的效率      #a_IB:鸟类捕食昆虫的效率
    #d_B:鸟类的自然死亡率              #e_BtI:蝙蝠从昆虫获得能量的效率
    # m_IBt:蝙蝠捕食昆虫的效率         #d_Bt:蝙蝠的自然死亡率

    dCdt = r_C * C - a_CI * I * C - d_C * C#农作物增长
    dWdt = r_W * W - h * W - d_W * W#杂草增长
    dIdt = e_IC * a_IC * C * I - m_IB * B * I - m_IBt * Bt * I - p * I#昆虫增长
    dBdt = e_BI * a_IB * I * B - d_B * B#鸟类增长
    dBtdt = e_BtI * m_IBt * I * Bt - d_Bt * Bt#蝙蝠增长

    return [dCdt, dWdt, dIdt, dBdt, dBtdt]


# 初始种群数量
C0 = 500  # 农作物
W0 = 500  # 杂草
I0 = 100  # 昆虫
B0 = 6  # 鸟类
Bt0 = 5  # 蝙蝠
y0 = [C0, W0, I0, B0, Bt0]

# 时间点
t = np.linspace(0, 200, 1000)

# 参数
'''r_C, a_CI, d_C, r_W, h, d_W, e_IC, a_IC, m_IB, m_IBt, p, e_BI, a_IB, d_B, e_BtI, d_Bt = params

    #r_C:农作物的自然增长率            #a_CI:昆虫对农作物的攻击率
    #d_C:农作物自然死亡率               #r_W:杂草的自然增长率
    #h:除草剂的效果系数                #d_W:杂草自然死亡率
    #e_IC:昆虫从农作物获得能量的效率    #a_IC:昆虫捕食农作物的效率
    #m_IB:鸟类捕食昆虫的效率           # m_IBt:蝙蝠捕食昆虫的效率
    #p:杀虫剂的效果系数                #e_BI:鸟类从昆虫获得能量的效率
    #a_IB:鸟类捕食昆虫的效率           #d_B:鸟类的自然死亡率
    #e_BtI:蝙蝠从昆虫获得能量的效率     #d_Bt:蝙蝠的自然死亡率'''


params = ()

# 求解微分方程组
solution = odeint(model, y0, t, args=(params,))

# 绘制结果
plt.figure(figsize=(10, 8))
plt.plot(t, solution[:, 0], label='农作物数量 (C)')
plt.plot(t, solution[:, 1], label='杂草数量 (W)')
plt.plot(t, solution[:, 2], label='昆虫数量 (I)')
plt.plot(t, solution[:, 3], label='鸟类数量 (B)')
plt.plot(t, solution[:, 4], label='蝙蝠数量 (Bt)')

plt.legend()

plt.xlabel('时间')
plt.ylabel('种群数量')
plt.title('生态系统种群数量随时间变化')
plt.show()